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High precision single-cluster Monte Carlo measurement of the critical 
exponents of the classical 3D Heisenberg model * 

C. Holm'^and W. Janke*^ 
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Arnimallee 14, D-1000 Berlin 33, Germany 

'^Hoclistleistungsrechenzcntrum, Forschungszentrum Jiilich 
Postfach 1913, D-5170 Jiilich, Germany 

We report measurements of the critical exponents of the classical three-dimensional Heisenberg model on simple 
cubic lattices of size with L — 12, 16, 20, 24, 32, 40, and 48. The data was obtained from a few long single- 
cluster Monte Carlo simulations near the phase transition. We compute high precision estimates of the critical 
coupling Kc, Binder's parameter U* , and the critical exponents v,(5/v,rj, and a/v, using extensively histogram 
reweighting and optimization techniques that allow us to keep control over the statistical errors. Measurements of 
the autocorrelation time show the expected reduction of critical slowing down at the phase transition as compared 
to local update algorithms. This allows simulations on significantly larger lattices than in previous studies and 
consequently a better control over systematic errors in finite-size scaling analyses. 
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1. INTRODUCTION 

The three-dimensional (3D) classical Heisen- 
berg model is one of the simplest spin models, 
and its critical behavior has been investigated 
by a variety of approaches. Recently Peczak, 
Ferrenberg, and Landau ^ (PFL) have under- 
taken a high statistics Monte Carlo (MC) study 
of this model on cubic lattices of sizes up to 
V = Li^ = 24"^, using standard Metropolis and 
multi-histogram techniques ||]. They found a 
value of Kc = 0.6929(1) that was in disagreement 
with previous estimates of the critical coupling 
Kc = J/ksTc- These estimates were derived from 
high-temperature series (HTS) expansion analy- 
ses D based on the Fade (Kc = 0.6924(2)) and 
ratio {Kc — 0.6925(1)) method, respectively, and 
more recent transfer-matrix (TM) MC investiga- 
tions D {Kc = 0.6922(2) and Kc = 0.6925(3)). 
The critical coupling is a non-universal parame- 
ter and from this point of view of no particular 
interest. Most estimates of universal critical ex- 
ponents, however, are biased and usually depend 
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quite strongly on the precise value of Kc- To clar- 
ify the above discrepancy we performed an inde- 
pendent high precision single-cluster MC study 
on large lattices of sizes up to 48^ . We could con- 
firm PFL's value of Kc and obtained measure- 
ments of the critical exponents that are in their 
accuracy comparable to the best estimates com- 
ing from field theory. 



2. MODEL AND ALGORITHM 

The Metropolis (pseudo) dynamics suffers from 
the severe problem of critical slowing down, which 
is the reason why we used as update algorithm 
the cluster algorithm in its single-cluster version 
H. One update in the single-cluster variant con- 
sists of chosing a random mirror plane and a ran- 
dom starting site, which is the germ of a cluster 
whose members are selected from adjacent sites 
by a Metropolis-like accept /reject criterion. From 
studies of related spin models it is known ||^ that 
this variant of the cluster algorithm is extremely 
efficient in three dimensions. 

The classical Heisenberg model is defined by 
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the Hamiltonian 

7i = J ^ [1 - • sj] (1) 

where J is the ferromagnetic coupling (J > 0). 
The sum runs over all nearest neighbour pairs 
(i, j), and sare three-dimensional unit spins that 
live on the sites i of a simple cubic lattice with pe- 
riodic boundary conditions. The continuous en- 
ergy range < i? < 3V was discretized into 90000 
bins, which is fine enough to ensure no significant 
discretization errors. We stored histograms in- 
stead of the whole data time series in order to 
save storage space. 

Our simulations were organized as follows. 
First, we did one run for each lattice size at 
Kq = 0.6929, the estimate of by PFL, 
and recorded the energy histogram Pko{E) 
and the microcanonical averages {{m'')){E) = 
Y.MPKo{E.M)m''/PK,{E), k = 1,2,4, where 
m = M /V = \fh\ is the magnitude of the magne- 
tization rfi = y ^(^) of a single spin config- 
uration. The temperature independent averages 
((m*^))(i?) can be computed by accumulating the 
values of m'^ in lists indexed by the associated 
energy bin of the configuration and normalizing 
at the end by the total number of entries in each 
bin, making it thus unnecessary to store the two- 
dimensional histogram Pkq{E,M). 

From the data of the run we could com- 
pute the approximate positions > Kq and 
< Kq of the (connected) susceptibility and 
the specific-heat peak maximum by reweighting 
techniques ||^. We then performed two more 
runs at K^ and K_ , respectively, again recording 
Pk{E) and {{m^)){E). This choice of the simula- 
tion points has the advantage that one automat- 
ically stays in the critical region since both K+ 
and K^ scale with L^^/^ , where v is the corre- 
lation length exponent. From this data we then 
computed three estimates 0^^\k), n — — ,0, + 
for all thermodynamic observables Ol of interest, 
and for any K value in the vicinity of K- , Kq, K+ 
by reweighting. The reweighting range was deter- 
mined by the energy value at which the energy 
histogram had decreased to a third of its maxi- 
mum. This ensured a high enough statistics of the 
histogram to allow the reweighting scheme to pro- 



duce reliable results. Furthermore we used block- 
ing to compute jackknife errors AO^"-* on C^"^ 
To get an optimized average of these three values 
that minimizes the relative error of the combined 
Ol{K) for each observable separately, we added 
the O^"^ with relative weight 1/(A0["^)2. All our 
runs contain at least 10000 x r measurements, 
where r is the integrated autocorrelation time 
of the susceptibility. The measured value of r 
turned out to be almost independent of the lattice 
size and to be very small, r < 2, in units of lat- 
tice sweeps that allow direct comparison with the 
Metropolis algorithm. This is of course expected 
for the single-cluster update. For the 48^ lattice 
we obtain a value of r which is about three orders 
of magnitude smaller than for the Metropolis al- 
gorithm. This explains why we could study much 
larger lattice sizes than PFL, and could still af- 
ford to have about ten times better statistics. 

3. RESULTS 

To determine Kc we first concentrated on 
Binder's parameter 

Ul{K) = 1-,^. (2) 

Asymptotically for large L, all curves Ul(K) 
should cross in the unique point {Kc,U*). The 
locations of the crossing points of two different 
curves Ul{K) and Ul'{K) depend on the scale 
factor h = L/L' , due to residual corrections to 
finite-size scaling (FSS). From the crossings with 
the L = 12 and L = 16 curves we obtained six 
and five data points, respectively, whose linear 
fits shown in Fig. 1 gave us the result 

A'c = 0.6930 ±0.0001, (3) 

and 

[/*= 0.6217 ±0.0008. (4) 

The critical coupling is found in excellent agree- 
ment with the value quoted by PFL, and also U* 
agrees very well with their estimate of 0.622(1). 
For comparison, a field theoretic e-expansion pre- 
dicts a 4% lower value of 0.59684... ||. 

To obtain an estimate for the correlation length 
exponent v we use that at K^ the derivatives 
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dUL/dK should scale asymptotically with L^/^. 
To calculate dUL / dK from our data, we took the 
thermodynamic derivative 



dUL 
dK 



{l-U){{E) 



{m?) (to*) 



(5) 



which is less vulnerable to systematic errors than 
a finite difference approximation scheme. In a 
log-log plot we find a perfect straight line fit 
(with goodness-of-fit parameter Q — 0.61 at Kc = 
0.6930) that yields 



0.704 ±0.006. 



(6) 



This value is in good agreement with PFL's mea- 
surement oi ly ~ 0.706(9) (determined by the 
same method, but a,t K — 0.6929), and with 
the field theoretic estimates derived from the re- 
summed (7-expansion |^ or from the resummed 
e-expansion see Table 1. The high quality of 
this fit (as well as of all other fits described below) 
shows that the asymptotic scaling formula works 
down to our smallest lattice size L = 12, indicat- 
ing that there is no need for confluent correction 
terms. 

Having now measured v we can get two more 
estimates for the critical coupling by assuming 
the FSS relation T^ax ^ + aL-^^" + ... for 
the location of the maxima of the specific heat 
C = V-^K^ ((^2) - (Ef), and the connected 
susceptibility x'^ = V K (^{m'^) — (m)^) . Using 
our value of = 0.704 we obtain from the linear 
fits shown in Fig. 1 the estimates Kc = 0.6925(9) 
(from Tc.„,, with Q = 0.80) and K^ = 0.6930(3) 
(from Ty,^^_^^ with Q ~ 1.0), respectively. These 
values are consistent with the crossing value (^, 
but have larger statistical errors. 

The ratio of exponents /3/i^ follows from the 
scaling of the magnetization, (to) oc L"^/''. From 
the linear least-square fit in a log-log plot of (to) 
versus L at Kr — 0.6930 we obtain the estimate 



/3/i/ = 0.514 ±0.001, 



(7) 



(with Q — 0.68) that is slightly lower than the 
value given by PFL, P/v = 0.516(3) (determined 
at K = 0.6929). To test by how much our result 
is biased by the value of Kc we have redone our 
analysis at X = 0.6929 [K = 0.6931). Here we 



Figure 1. Estimates of the critical temperature 
Tc, coming from Binder's crossing method and 
from the scaling of Tc^^^ and T^c 



get the slightly higher (lower) value of 0.519(1) 
(0.509(1)). The quality of the fits, however, is 
much worse, namely Q = 0.30 {Q = 0.31), which 
we interpret as support for our estimate of Kc. 
We rely on the goodness-of-fit parameter since 
visually it is impossible to make a distinction be- 
tween these fits when plotted on a natural scale. 
It should be emphasized that even these slight 
variations in the estimate of the critical coupling 
significantly change the estimate of the exponent 
ratio in a way that clearly dominates the statis- 
tical errors, making it necessary to have an accu- 
rate estimate of Kc. 

Relying on the scaling law rj = 2 — j/i^ =- 
2(3 /v — 1 we estimate 



77 0.028 ±0.002. 



(8) 



One can get independent estimates of 77 by direct 
measurements of x = VK{m'^) and x'^ our best 
estimate of Kc = 0.6930, where they both should 
scale like L^~^. From linear fits we obtained rj = 
0.0271(17) (Q = 0.78) and 77 = 0.0156(44) (Q 
0.69), respectively. Finally, analyzing the FSS 
behavior of the susceptibility maximum, Xma->^ 
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Table 1 

Selected sources for K^^ U*, and the critical exponents of the classical 3D Hciscnbcrg model. 



method 



U* 



V 



ajv 



^r-expansion 
e-expansion |7|J 
HTS [P 
MC 

MC (this work) 




0.6929(1) 
0.6929(1) 
0.6930(1) 



0.59684 

0.622(1) 
0.6217(8) 



0.705(3) 

0.710(7) 

0.712(10) 

0.706(9) 

0.704(6) 



0.517(6) 

0.518(11) 

0.513(50) 

0.516(3) 

0.514(1) 



0.033(4) 

0.040(3) 

0.034(42) 

0.031(7) 

0.028(2) 



-0.163(12) 
-0.183(28) 
-0.191(40) 
-0.167(36) 
-0.159(24) 



l?--^, we estimate = 0.0231(61) (Q = 0.60). 
Notice that all MC estimates are lower than the 
field theory values, which are collected in Table 1. 

Similarly, using the hyper-scaling law ajv — 
2/v — ?i we obtain for the specific-heat exponent 

a/i^ = -0.159 ±0.024. (9) 

Measurements of the specific heat are difficult to 
analyze directly, because a is negative, which im- 
plies a finite, cusp-like singularity. We tried a 
three-parameter fit of the form Cmax = a — bL"/'^. 
The result a/v = -0.33(22) (Q = 0.69) is com- 
patible with eq. (H), but has large statistical er- 
rors. Another way of testing eq. (||) is to assume 
the predicted value of a/v and to fit only the pa- 
rameters a and b. The resulting fit turned out to 
be of almost equally good quality |l^ . 

To summarize, using the single-cluster MC up- 
date for the classical 3D Heisenberg model on sim- 
ple cubic lattices of size up to 48'^, we obtained 
high-precision data. Using multi-histogram tech- 
niques we optimally combined the data and per- 
formed a fairly detailed FSS analysis. Qualita- 
tively, our main result is that the asymptotic FSS 
region sets in for small lattices sizes, L « 12. 
Quantitatively, our value for the critical coupling, 
Kc = 0.6930(1), is significantly higher than es- 
timates from old HTS expansion analyses and 
TM MC methods, but is in almost perfect agree- 
ment with the MC estimate reported recently by 
PPL [P, with new analyses of longer HTS expan- 
sions |ll[ , and with recent MC simulations in the 
high-temperature phase ||l^. Our results for the 
two basic critical exponents, ly = 0.704(6) and 
/3 = 0.362(4), are in good agreement with field 
theoretic predictions. Direct measurements of rj 
give good agreement with the scaling prediction 
when the scaling of x at Kc is considered. Using 



at Kc or Xmax: however, the situation is less 
clear. In the case of a/v, its negative value causes 
numerical problems, which result in large statis- 
tical errors. Table 1 lists our measured values 
and their scaling implications for rj and a/v. For 
comparison, various other sources for the critical 
exponents are added. More details of this study 
can be found in refs. |p^ , |l2| . 
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